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ABSTRACT 

A recently reported solution for stationary stability of a thermo- 
solutal system with Soret diffusion is re-derived and examined using a 
symbolic computational package. Symbolic computational languages 
are well suited for such an analysis and facilitate a pragmatic approach 
that is adaptable to similar problems. Linearization of the equations, 
normal mode analysis, and extraction of the final solution are performed 
in a Mathematica ® notebook format. An exact solution is obtained for 
stationary stability in the limit of zero gravity. A closed form expres- 
sion is also obtained for the location of asymptotes in relevant param- 
eter, (Sm c ,Ma c ), space. The stationary stability behavior is conveniently 
examined within the symbolic language environment. An abbreviated 
version of the Mathematica ® notebook is given in the Appendix. f 
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surface Sherwood number, 
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Soret coefficient, 


temperature 
time 
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surface tension 
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u disturbance velocity vector 

V y-component velocity 

W z-component velocity 

w z-component disturbance velocity 

x,y,z Cartesian coordinates (see Fig. 1) 


Greek Symbols 


NOMENCLATURE 

C concentration 
c disturbance concentration 
D m n diffusivity elements surface tension 
Dp Yi 

Dm Dufour coefficient, — 

Dll 72 

d depth of fluid layer 
h surface heat transfer coefficient 
h s surface mass transfer coefficient 
J mass flux from surface 

YidAT 

Ma Marangoni number, — 

DuM 

hd 

Nu surface Nusselt number, 

P c p D ii 

Pr Prandtl number, v/D u 
AP pressure difference 
Ap disturbance pressure difference 
Q heat flux from surface 


a wavenumber 

y surface tension variation 

9 disturbance temperature 

4> normal mode temperature 
X eigenvalue 

p dynamic viscosity 

v kinematic viscosity 

p density 

o surface tension 

x diffusivity ratio, D 22 /Dj j 

X normal mode concentration 

Subscripts 

0 value at z = 0 

1 value at z = 1 

b basic state 

c critical value 

oo asymptotic value 


tThe unabridged Mathematica ® notebook is available from the author or on 
the internet at http://zeta.lerc.nasa.gov/6712/people/skarda.htm 
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INTRODUCTION 

This paper is intended to achieve three goals. First the paper and 
an associated Mathematical notebook document the derivation details 
of an exact solution for stationary stability of the Marangoni-Benard 
system with Soret diffusion (Skarda et. al., 1998). Second, it is hoped 
that the problem described herein might serve as a “hands on learning 
toor to augment teaching or understanding of linear stability. The indi- 
vidual steps of a linear stability analysis and their relationship to each 
other can effectively be examined within a symbolic computational 
environment. Finally, several features of symbolic computational 
languages such as Mathematical (Wolfram, 1996), are demonstrated 
within the notebook that can greatly benefit similar analyses. In fact, 
various elements of the notebook are directly, or with minor modifica- 
tion, applicable to a wide range of stability problems. For example the 
brief linearization function we construct is quite general and can be 
applied to linearize systems of nonlinear partial differential equations 
(p.d.e.) of arbitrary order and degree, as well as the associated bound- 
ary conditions. Similarly, the normal mode substitutions apply to most 
autonomous systems, where the boundary conditions can also be satis- 
fied by the normal mode solutions. 

The physical problem under consideration is that of binary fluid 
convection with soret diffusion. Hurle and Jakeman (1969 and 1971) 
showed that these binary systems can be quite sensitive to the Soret 
effect resulting in either undesirable flows due to trace amounts of im- 
purities or desirable mixing by providing a small concentration of 
solute. Binary fluid systems are also instructive because they can be 
represented as relatively simple fluid-dynamical systems with constant 
coefficients (Jacqmin, 1996) and zero base flow. The convection dis- 
plays an abundance of interesting phenomena as demonstrated by sev- 
eral theoretical and experimental investigations, e g., Kolodner et. al. 
(1989), Bensimon et. al. (1990), Glazier and Kolodner (1991), Cross 
(1986), Knobloch and Moore (1988). In these studies, flow is driven 
due to density variations. More recent investigations have examined 
the Soret effect in systems where surface tension gradients at the free 
surface drive flow (Chen and Chen 1994, Castillo and Velarde, 1978, 
and Bergeon, et. al., 1994, Skarda, et. al., 1998). Such phenomena is 
important to proposed containerless processing applications in 
microgravity. 

We discuss the linear instability of double diffusive convection 
with the Soret effect that is driven by surface tension variation along 
the free surface. An imposed temperature difference across the layer, 
AT b , induces a concentration difference, AC b . We rescale the problem 
in the absence of buoyancy which leads to a concise representation of 
neutral stability results in and near the limit of zero gravity. An exact 
solution is obtained for stationary stability. One important consequence 
of the exact solution is the validation of published numerical results in 
the limit of zero gravity. Moreover, the precise location of asymptotes 
in relevant parameter, (Sm c ,Ma c ), space are computed from exact solu- 
tions. Stability behavior is briefly examined using the exact solution. 
The process of linearizing the equations, deriving the normal mode 
equations, and obtaining an exact solution for neutral stability are illus- 
trated in the Appendix in Mathematical notebook format. 


DEVELOPMENT OF EQUATIONS 

We consider the unbounded cross-doubly diffusive fluid layer of 
depth d as shown in Fig. 1. Buoyancy is neglected, and onset of con- 
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Figure 1 . — Unbounded double diffusive fluid layer. 


vection due to surface tension variation is examined. The basic or un- 
perturbed state is one of no flow, U b = 0, an imposed temperature pro- 

file across the layer, T b (z) = T b (0)- AT b — , and an induced 

d 

concentration gradient across the layer, AC b =-^^-AT b . Difference 

D 22 

quantities of the form, Ay b , are defined as Ay b = y b (0)-y b (d). 

The governing equations we wish to linearize are: 


vu = o 

(1) 

P-5^ = -VP + hV 2 U 
Dt 

(2) 

DT - o 

— = D u V 2 T + D 12 V 2 C 

(3) 

— =d 21 v 2 t+d 22 v 2 c 

Dt 21 22 

(4) 


The dynamic viscosity, ji, and diffusivity elements, D mn , are assumed 
constants. The boundary conditions are: 
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where V tI = i — + j — . Surface tension, a, is approximated as a lin- 
3x dy 

earized function of T and C, a = a 0 - Yj(T - T b ) - y 2 (C - C b ), where y x 
and y 2 are the temperature and concentration variations with sur- 
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cells for the linearized equations after evaluating the solution forms, 
Eqs. (A.9a) to (A.9c). The resulting normal mode equations are also 
given in (Skarda et. al., 1998). The energy and species equations are 
coupled by their cross diffusive terms which also occur in the flux bound- 
ary conditions at z = 1 . The momentum equation is coupled with both 
energy and species equations through the tangential stress condition, 
Eq. (A. 8b). 

The eigenvalue, X in Eqs. (9a) to (c) is complex. If the real part of 
X is positive, a disturbance in the system grows, and if the real part is 
negative, the disturbance decays. If X = 0, the disturbance persists un- 
changed in time and the system is in a state of neutral stability. By 
setting X = 0, we can decouple the energy and species equations and 
directly solve the following normal mode equations for stationary 
stability. 


J = j(T b ,C b ) + |t (C(l)-C,) (7b) 

Substituting Eq. (7a) into Eq. (6d), and Eq. (7b) into Eq. (6e) for Q and 
J, respectively, facilitates linearization of the flux conditions. 

The dependent variables are written in terms of the base flow and 
perturbation variables, U = u, T = T b + 0, C = C b + c, and AP = AP b + 
Ap. Substituting these into the governing equations and boundary con- 
ditions, and retaining only the first order terms leads to the linearized 
system. Exploiting the pattern matching and limit evaluation capabili- 
ties of Mathematical, a succinct function was constructed to linearize 
a nonlinear p.d.e of arbitrary order and degree. This linearization func- 
tion, Linearize [equation], is given in the Appendix. Our problem was 
easily linearized, symbolically, as is illustrated within the notebook in 
the Appendix. 

As is common practice (Chandrasekhar, 1981) the curl operator is 
applied twice to the momentum equation. This operation could also be 
performed using the Mathematical vector calculus operators. 

Following the classical treatment for the Benard type problems, 
the linearized equations are now nondimensionalized (Chandrasekhar, 
1981). Here, they are nondimensionalized using the reference values of 
Skarda et. al. (1998) for length, velocity, time, temperature, and con- 
centration which are defined as 


0 = Pr(D 2 -a 2 ) 2 w 

(8) 

(d 2 -a 2 j<() = d < j ) Ma w 

(9) 

(D 2 -a 2 )x = d x Ma w 

(10) 

where 


, T+SmDm , , 1+t 

d A = and d 7 = 

v x — xSmDm A x — xSmDm 


The corresponding normal mode boundary conditions are given 
below: 

At the lower surface, z = 0, 

w(0) = 0 

(11a) 

Dw(0) = 0 

(lib) 


d. 


, — , AT b , and - AT b . 

D u D 


4>(0)=o 


(lie) 


The resulting nondimensional parameters are the Prandtl number, Pr, 
diffusivity ratio, x, Dufour coefficient, Dm, Soret coefficient, Sm, 
Marangoni number, Ma, surface Nusselt Number, Nu, and surface 
Sherwood number, Sh as defined in the nomenclature. The complete 
set of dimensionless disturbance equations are given by Eqs. (A.4) to 
(A.8d) in the Appendix. 

A normal mode analysis is performed by assuming solutions of 
the form (u, 0, c) = (w(z), (|)(z), ^(z))e^ l+l(a * x+a y y ). In Mathematical 
this is conveniently done by defining the solutions as the functions, 
Eqs. (A. 9a) to (A.9c) in the Appendix. 

These solutions are substituted into the set of linearized Eqs. (A.4) 
to ( A.8d) leading to a 4 th order and two 2 nd order ODE’s. The normal 
mode equations are easily obtained by evaluating the Mathematical 


<K0) = 0 (11c) 

At the upper surface, z = 1, 


w = 0 

(12a) 

-D 2 w = a 2 (<(> - Smx) 

(12b) 

D6 - Sm Dm D/ + Nu 0 = 0 

(12c) 

-D<(> + Dx + Sh x = 0 

(1 2d) 
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STATIONARY SOLUTION 

The solution forms of w, and x are easily determined from in- 
spection of Eqs. (8) to (10) and are given in the Appendix by Eqs. (A. 12) 
to (A. 14). Determination of the integration constants for this problem 
is tedious. Evaluation of bl and cl from Eqs. (A.21) and (A.22) should 
confirm this claim. A highly desirable feature of symbolic languages is 
their ability to minimize much of the computational drudgery, verify 
results obtained by other means, and often simplify the final expression 
to a tractable form. The use of Mathematical to solve for integration 
coefficients, bl and cl, certainly minimizes drudgery, and we will find 
below that the final form of the neutral stability solution is simplified 
by Mathematical to a tractable form. An expression for neutral stabil- 
ity is obtained by applying the tangential stress balance, Eq. (12b). This 
yields the following exact solution, Eq. (13), for stationary stability of 
the double diffusive Soret problem (Skarda et. ai., 1998). For general- 
ity the Dufour diffusivity term is also retained in the energy equation. 


Maj 

( Sm 

r, 

1 T ) 

\ 

-1 

) 

a cosh a + | 

(v sm " sh ] 

sinhaj 

l 

2i 

^(Sm Dm - 

- l)a 2 (l + cosh 2a) + Nu Sh(l - cosh 2a) -(Nu + Sh)asinh2a)J 


8a 2 - 4asinh2a 

+ ^ = 0 (13) 

sinh3a- 3sinha-4a cosh a 


RESULTS 

The neutral stability curves shown in Figs. 2 and 3 are computed 
directly from Eq. (13). The abscissa and ordinate are the wavenumber 
and the Marangoni number, respectively, and the eigenvalue, X, is zero 
along each curve. For Ma values above a neutral stability curve (corre- 
sponding to given Sm and T values), infinitesimal disturbances grow in 
time and the fluid layer is said to be unstable. The layer is linearly 
stable for Ma values below the appropriate neutral stability curve. In 
Fig. 2, the curves are shifted upward for increasing values of Sm thus 
increasing Sm has a “stabilizing effect.” In Fig. 3, decreasing T values 
are destabilizing, shifting the neutral stability curves downward. We 
emphasize that these are stationary stability results; however, oscilla- 
tory stability will occur before stationary stability (Ma values below 
the stationary values) for certain ranges of Sm values (Chen and Chen, 
1994, Skarda et. al., 1998). 

The minimum value of Ma on each curve in Figs. 2 and 3 is called 
the critical Marangoni number, Ma c , and is of interest since onset of 
instability typically occurs at this value. In Fig. 4, stationary stability 
boundaries, the loci of points consisting of Ma c values, are shown in 
(Sm c ,Ma c ) space for different values of T. The overall behavior is quite 
similar to the Rayleigh-Benard stationary stability results presented in 
(Sr c ,Ra c ) space by Hurle and Jakeman (1971) where the Soret coeffi- 
cient, Sr, and Rayleigh number, Ra, are defined in the usual manner 
(Hurle and Jakeman 1971). Asymptotic behavior occurs at a finite Sm c 
value and Ma c is driven to zero as ISrnJ — > «>. The exact location of the 
asymptote observed in Fig. 4 is also determined from our solution, 
Eq. (13), and is given in Skarda et. al., (1998). For flux boundary con- 
ditions often applied to the Soret problem, Nu = Sh = 0, and Dm = 0 the 
location of the^ asymptote reduces to a simple function of T, 

Sm^ = ^1 + — j . This result was found to be identical to the asymp- 



Figure 2. — Neutral stability curves, each curve for constant 
Sm values, t = .05, Dm = 0, Nu = 0, Sh = 0. 



Figure 3. — Neutral stability curves, each curve for constant 
t values, Sm = -.02, Dm = 0, Nu = 0, Sh = 0. 


tote location expression obtained by Hurle and Jakeman (1971), but in 
(Sr c ,Ra c ) space, for the case of the buoyancy induced Soret problem 
with free-free surfaces (Skarda et. al. 1998). As before oscillatory 
instability which must be computed numerically, is not considered in 
this analysis, but is discussed in the previously mentioned references. 


CONCLUSIONS 

The exact solution for stationary onset of convection for 
Marangoni-Benard instability with Soret effect was successfully 
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Sm c 


Figure 4. — Stationary stability boundaries, each curve for 
constant t values, a c = 1 .99, Dm = 0, Nu = 0, Sh = 0. 


derived using a symbolic computational package, Mathematical . The 
results agree with the “hand-derived” form of the solution recently 
reported in the literature. The associated notebook permits users to 
reconstruct the steps leading to the solution. The exact solution is also 
provided in a form to allow immediate examination of neutral stability 
behavior and exploration of various parametric effects. While exact 
solutions are uncommon, various parts of the notebook can be applied 
with little modification to assist with aspects of other linear stability 
problems. 
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APPENDIX 


Exact Solution For Stationary Stability Double Diffusive Marangoni Problem With 
Soret Diffusion 

Cross-Diffusion Terms Are Included in Governing Equations & in the Flux BC’s at z=1 
A Conductive/Permeable Base Is Considered at z=0 


<<Calculus 9 VectorAnalysis ' ; 

SetCoordinates [ Cartesian [x, y, z] ] ; 

(* Define linearization function *) 

Linearize [equation^] := Limit [Expand [ (Expand [equation] /. e n - -> 0) / e] , e -> Infinity] 


■ Governing Equations - Momentum, Energy, Species 

UMomEq = (* A . 1 *) 

tf[x # Y * t] + (U[x, y, z, t], V[x, y, z # t], W[x # y, z, t]} . Grad[U[x, y, z, t] ] - 

d x P[x, y, z, t] + Laplacian [ U[x, y, z, t] ]; 

EnergyEq = (* A .2 *) 

T[x, y, z, t] + {U[x, y, z, t] , V[x, y, z, t] , W[x, y, z, t]} .Grad[T[x, y, z, t] ] - 

Du Laplacian [ T [x, y, z, t] ] - D 12 Laplacian[ S [x # y, z, t] ]; 

SpeciesEq = (★ A. 3 *) 

d t S [x, y, z, t] + (U[x, y, z, t] , V[x, y, z, t] , W[x, y f z, t] } . Grad[ S [x, y, z, t] ] - 
D 2 i Laplacian [ T [x, y, z, t] ] - D 22 Laplacian [ S [x, y, z, t] ] ; 

U[x_, y_, z_, t_] : = 0+ e u [x, y, z, t] ; 

V[x_, y__, z_, t_] : = 0 + 6 v[x, y, z, t] ; 

W[x_, y_, z_, t_] := 0+ ew[x, y, z, t] ; 

P[x_, y_, z_, t_] := 0 + ep [x, y, z , t j ; 

T[x__, y_, z_, t_] := T b [x, y, z, t] + e6[x, y, z, t] ; 

S[x_, y_, z_, t_] : = S b [x, y, z, t] + ex[x, y, z, t] ; 

Linearize [UMomEq] ; 

Linearize [EnergyEq] ; 

Linearize [SpeciesEq] ; 


■ Linearized Disturbance Equations 


MDisturbEq = d t Laplacian [ u [x, y, z, t] ] - Pr Biharmonic [u [x, y, z, t]] (* A .4 *) 

EDisturbEq = d t 6[x, y, z, t] - Laplacian [6 [x, y, z, t]] + (★ A. 5 *) 

Sm Dm Laplacian [c [x, y, z, t] ] - Mau[x, y, z, t] 

SDieturbEq = d t c [x f y, z, t] + z Laplacian [0[x, y, z, t] ] - (* A.6 *) 

z Laplacian [c [x, y, z, t] ] - Mau[x, y, z, t] 

■ Boundary Conditions 


■ at z=0 


u[x, y, 0 # t] ; 
d%.z u[x, y, 0, t] ; 
0[x, y , 0 # t] ; 
c [x, y, 0, t] ? 


(★ A . 7 , a, b, c, d *) 
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■ at z=l 


u[x, y, 1, t] ; 

d*,z U[x, y, 1, t] + a XfX u[x, y, 1 # t] + d y , y u[x # y, 1, t] ; 

d z ©[x, y# 1# t] + Dm 6 Z c [x, y, 1, t] +Nu0[x, y, 1, t] ; (* A, 8 a,b,c,d *) 

Sm a z 0[x, y, 1 , t] + r d z c [x, y, 1, t] + Sh c [x, y, 1, t] ; 


■ Assume Normal Mode Solutions of The Form 

u[x_, y_, z_, t_] : = w[z] Exp[At + I (a x x + a y y) ] ; 

e [x / y_, z_, t_] : = <t>[ z] Exp [At + I (a x x + a y y) ] ; (* A . 9 a,b, c *) 

c[x_, y_, z_, t_] := x[z] Exp [At + I (a x x + a y y) ] ; 


■ Normal Mode Equations 

(D 2 - or 2 ) 2 w[z] = 0 

(D 2 - o 2 )^] = -d# Ma w[z] 

(D 2 - a 1 )xlz] = -d x Ma w[z] 

(* where *) 

d<p = (t + Sm Dm) / (i-tSm Dm) ; dx = (l+r)/(r-r Sm Dm) 


(* A.10 *) 
(* A. 11 *) 
(* A.12 *) 


■ Solutions For velocity, Temperature, and Concentration Are Expressed As: 


w[x__] := al Sinh[a x] + a2 Cosh[a x] + (* A. 13 *) 

a3 x Sinh [a x] + a4 x Cosh [a x] 

^[x_] := bl Sinh[a x] + b2 Cosh [a x] + 

b3 x Sinh[a x] + b4 x Cosh [a x] + (* A. 14 *) 

b5 x A 2 Sinh [a x] + b6 x A 2 Cosh[a x] 


X[x_] := cl Sinh[a x] + c2 Coshfa x] + 

c3 x Sinh[a x] + c4 x Cosh[a x] + (* A. 15 *) 

c5 x A 2 Sinh[a x] + c6 x A 2 Cosh[a x] 


■ Solve for velocity coefficients, a2, a3, and a4 in terms of al 


Boundary Condition, w[0] = 0, w’[0] , and w[1]=0 Give 


a2 = 0; a3 = al (a (Cosh [a] /Sinh [a] ) -1) ; a4 = -al a; (* A.16a,b,c) 


■ Solve for particular temperature & species solutions 

Solve for Particular Solution of <j> (temperature), 
ie solve for b3, b4, b5, b6 (also c3,c4,c5,and c6 values) 


D[b3 al z Sinh [a z] + b4 al z Cosh [a z],{z,2>] - 

a A 2 (b3 al z Sinh[a z] + b4 al z Cosh[a z]) + (* A. 17 *) 

D[b5 al z A 2 Sinh[a z] + b6 al z A 2 Cosh[a z],{z,2}] - 

a A 2 (b5 al z A 2 Sinh[a z] + b6 al z A 2 Cosh[a z] ) /. al->l /. d^>->l 
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Values of b3, b4, b5, b6 (or c3, c4, c5, or c6) Determined by comparison 
with Cosh, Sinh, etc Terms of ODE le w[z] Terms 


b3 

* al d 0 Ma 

(a 

Coth [a] -1) / (4 

a A 2) ; 

b4 = -3 al a<t> Ma 

/ 

(4 

a) ; 

<* 

A. 18a, b 

*) 

b5 

= al d0 Ma / 

4; 

b6 = -al d0 

Ma (a 

Coth [a] -1) / (4 a); 




(* 

A. 18c, d 

*) 

c3 

= al dx Ma 

(a 

Coth [a] -1) / (4 

a A 2) ? 

c4 = -3 al dx Ma 

/ 

(4 

a) ; 

(* 

A. 19a, b 

*> 

c5 

■ al dx Ma / 

4; 

c6 = -al dx 

Ma (a 

Coth [a] -1) / (4 a) ; 




(* 

A. 19c , d 

*) 


■ Solve For Homogenous Coefficients - Energy & Species Equations 


(* Bottom surface, z=0, is conductive and permeable, therefore: *) 

b2 = 0; c2 = 0; (* A.20a,b *) 

Now apply the flux conditions at the top, z=1 and solve for bl and cl 


♦fluxbc = *'[1] + Nu 0[1] + (* A. 21 *) 

-Sm Dm x’ [1] ; 

Xfluxbc * - 4>' [ 1] + (* A. 22 *) 

x' [1] + Sh X[1J ; 


Solve For bl & cl simultaneously from flux conditons 


bl = . ; cl = . ; 
blcl m Flatten [ 

{bl,cl> /. Solve [ {0fluxbc==*O ,xfluxbc==0} , {bl, cl} ] ]; (* A. 23 *) 

bl = blcl [[I]]; cl = blcl[[2]]; 


■ Substitute w, <f>, and ^ into tangential stress boundary condition and obtain 
stationary stability solution 


(* The remaining coefficient, al will cancel out in boundary condition so we’ll set it to 1 *) 
01 = 0[1] /. al->l;xl = XtH / • al->l; 
d2wl = Simplify [Expand [w' ' [1] ] /. al->l ]; 


Substitute above into Tangential stress condition at z=1 


tstress = a A 2 01 - Sm a A 2 xl + d2wl; ( * A. 24 *) 

Resolve <p Term of Tangential Stress Equation - Multiply by 4 Sinh[a] to be Consistent With Literature 


4coef = 4 Simplify [a A 2 01] Sinh [a] 

Ma (a (Dm Sm + r) Cosh [a] + Sh t Sinh [a] ) (a 3 Cosh [a] - Sinh [a] 3 ) 
r ( (-1 + Dm Sm) a 2 Cosh[a] 3 - (Nu + Sh) a Cosh[a] Sinh[a) - Nu Sh Sinh[a] 2 ) 

0Term = 0coef /0coef [ [5 J ] 

Ma (a (Dm Sm + c) Cosh[a] + Sh z Sinh[a] ) 
r ( (-1 + Dm Sm) a 2 Cosh[a] 2 - (Nu + Sh) a Cosh [a] Sinh [a] - Nu Sh Sinh[a] 2 ) 


(* A .25 *) 


(* A .26 *) 
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Resolve* Term of Tangential Stress Equation - Multiply by 4 Sinh[a] to be Consistent With Literature 


xcoef = 4 Simplify [a A 2 xH Sinh[a] 

Ma (a (1 + z) Cosh [a] + Nu Sin h[a]) (a 3 Cosh [a] - Sinh [a] 3 ) 

z ((-1 + Dm Sm) a 2 Cosh [a] 2 - (Nu + Sh) a Cosh [a] Sinh [a] - Nu Sh Sinh [a] 2 ) 

XTenn = -Sm xcoef /xcoef [[ 5 ] 3 

Ma Sm (a (1 + z) Cosh [a] + Nu Sinh [a] ) 
z ( (-1 + Dm Sm) a 2 Cosh [a] 2 - (Nu + Sh) a Cosh [a] Sinh [a] - Nu Sh Sinh [a] 2 ) 

d2wcoef = Expand [ 4 Sinh [a] d2wl ]/tf>coef[[5]] 

8 a 2 - 4 a Sinh [2 a] 
a 3 Cosh [a] - Sinh [a] 3 


■ Collect the <f> and^ Terms into one term containing the soret coefficient Sm, 
d2wcoef, generate Neutral Stability Solution. 

NeutralSoln = Simplify [ (Numerator [^Term] + Numerator [xTerm] ) ] / 
Denominator [0Term] + Numerator [d2wcoef ] /Denominator [d2wcoef ] 

Ma (a (Sm (-1 + Dm - r) + r) Cosh [a] + (-Nu Sm + Sh t) Sinh [a] ) 
z ( (-1 + Dm Sm) a 2 Cosh [a] 2 - (Nu + Sh) a Cosh [a] Sinh [a] - Nu Sh Sinh [a] 2 ) 

8 a 2 - 4 a Sinh [2 a] 
a 3 Cosh [a] - Sinh [a] 3 


(* A .27 *) 
(★ A .28 *) 
(* A .29 *) 

Then with 

(* A. 30*) 
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